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Abstract 

We establish links between three different models for multifragmentation: 
the percolation model, the lattice gas model and the statistical multifrag- 
mentation model. There are remarkable similarities between the lattice gas 
model and the statistical multifragmentation model. For completeness, we 
also compare with a model based upon classical molecular dynamics which 
gives rather different results. 

PACS numbers: 21.60.-n, 21.10.-k 
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I. INTRODUCTION 



In this paper we establish relationships between three models commonly used for the 
theoretical treatment of nuclear multifragmentation. These are (1) the percolation model 
|Jl|J^, (2) the lattice gas model 0^ and (3) a class of models that we generically call the sta- 
tistical multifragmentation model (SMFM) A great deal of work in the last approach 
has been done. A complete description can be found in Ref. |^ and we will essentially follow 
the SMFM model presented there. Both the percolation and statistical multifragmentation 
model have been used for years. The lattice gas model is of more recent origin and can 
actually bridge the gap between those two seemingly totally different approaches. 

The connection between the percolation and the lattice gas model was already elaborated 
upon when the lattice gas model was introduced. It is included here for completeness. The 
main objective of the present paper is to point out that even though the lattice gas model 
and SMFM start off with two seemingly very different premises, there are many common 
features between the two models. Prompted by this similarity, we are driven to try and draw 
a. P — V diagram for the SMFM and to obtain the parameters of its critical behaviour, as 
for the lattice gas. We view this current work as mostly theoretical and we shall make only 
brief references to experimental results. For completeness we also make some comparisons 
with a classical molecular dynamics model that is in some use at the present time. 

II. THE LATTICE GAS AND PERCOLATION MODELS 

The lattice gas model has been described before ^j^. For completeness we describe some 
of its features that are relevant for this work. 

The lattice gas approach is viewed in our context as a modeling tool for the later stages 
of nucleus-nucleus collisions. For the case of n nucleons at this nuclear disassembly phase, 
the lattice gas theory consists of placing those n nucleons at lattice sites where N > n 
(when N = n the nucleus is at its normal volume). Each lattice site can be populated by 
one nucleon or none. In the lattice gas model we can not squeeze the nucleus to lower than 
normal volume. Later we will see that this feature is common to the SMFM also. The ratio 
N/n is equal to po/p/ where po is the normal nuclear density =0.16 fm~^ and pj is the 
freeze out density. The ratio N/ n is a parameter of the model and previously we have found 
that the value 2.54 provided a good fit with data 0. Two nucleons in adjacent lattice sites 
will interact with each other. We have previously taken this interaction to be attractive 
and equal to e, irrespective of whether it is between like particles (p-p and n-n) or between 
unlike particles (n-p interaction). We shall still adopt this argument for the time being but 
we shall later improve on this point. To generate an event we have to put nucleons in sites at 
a given temperature, taking into account that it is energetically preferable to put a nucleon 
adjacent to an occupied site rather than adjacent to an empty one. This Monte-Carlo 
sampling involves a Hamiltonian and is more involved than in a pure site percolation model. 
As can be readily guessed, the temperature is related to the excitation energy per nucleon. 
Once the nucleons have been placed in lattice sites we have to specify the momentum of 
each nucleon. The momentum of each nucleon is generated by a Monte-Carlo sampling of a 
Maxwell-Boltzmann distribution for the specified temperature. Thus at a given temperature 
we can generate an event in which the potential energy is obtained from the number of n — n 
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bonds and the total kinetic energy is simply the sum of kinetic energies of the n nucleons. 
We now recall how clusters are generated in the model. Once we have put in the nucleons in 
their lattice sites and generated their momenta it is easy to calculate the cluster distribution. 
Two nucleons in adjacent sites will form a cluster or be part of the same cluster if the relative 
kinetic energy of the two nucleons is insufficient to overcome the attraction: pf./2fi + e < 0. 
Here /i is the reduced mass, and e is a binding energy parameter from which the binding 
energy per nucleon can be calculated. Using this prescription, clusters in each simulation 
are easily found. It is also clear that in general, clusters do not have to appear in their 
ground state. They can be in an excited state. A comparison of yields of clusters calculated 
with this prescription was found in good agreement with many experiments Since the 

lattice gas model has a Hamiltonian, it is possible in principle to obtain an equation of state. 
The P — V diagram for the lattice gas model was drawn in ^ for an infinite system using the 
Bragg- Williams approximation. Results from improved Bethe-Peierls approximation were 
also shown. Clear indications of thermal critical phenomena are seen. It will be useful to 
remember here that the P — V diagram for lattice gas in the mean field approximation is 
similar to that of a Van der Waals gas. Improved calculations put the critical temperature 
at 1.1275|e|. The critical density is po/2 [p!0|] . 

The similarities of the lattice gas model with the percolation model are easy to see. For 
the latter we restrict ourselves to bond percolation model which is normally used for compar- 
ison with experiments. In terms of observables in heavy ion reactions, the bond percolation 
model is almost identical with the lattice gas model in the limit N = n. In the lattice gas 
model two nearest neighbours are bound to each other if + e < 0. This probabil- 

ity decreases with increasing temperature. In the bond percolation model the probability 
that bonds between nearest neighbours operate is given by a probability p which similarly 
decreases with increasing temperature. Indeed it was shown that the phenomenological 
parametrization for p used in bond percolation calculation makes the bonding probability 
same in both the models |]TT[]. For a certain value of p, a critical percolation cluster will 
appear. In a bond percolation model there is a concept of temperature but the volume is 
not a variable. Thus the concept of pressure is also absent. This is still a very valuable 
model since the fragment distribution, which is the experimental observable, is dictated by 
percolation characteristics. It is evident that the lattice gas model possesses the percolation 
model as a subset. In addition to containing all the percolation phenomena not only at nor- 
mal volume but also at larger volumes, it also has thermal phase transitions. The concepts 
of temperature, pressure and volume are all vahd concepts. An equation of state curve can 
be drawn for the lattice gas model but the percolation model is missing a parameter and we 
can not thus plot a traditional equation of state. 

III. THE STATISTICAL MULTIFRAGMENTATION MODEL 

Our arguments follow those presented in Ref. |^ except that since we are not as much 
interested in comparing with data as with connections between different models, we will 
make some simplifications. One assumes that equilibration has taken place in a freeze-out 
volume which is greater than normal nuclear volume. In this volume composites of different 
mass numbers appear. The volume of any of these composites will be its normal nuclear 
volume, i.e., fc/po where k is the mass number of the composite. It then follows that if the 
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total number of nucleons that are fragmenting is n, the excluded volume Vex is n/po and 
the available volume for thermalization is Vf = V — Vex- The concept of excluded volume is 
also present in the lattice gas model. Apart from the excluded volume the lattice gas model 
has a potential, as evidenced by the attraction e between nearest neighbour nucleons. This 
results in binding energy for each cluster. In the SMFM the specific two body interaction 
is not written down but there is clearly a potential energy and it makes its presence felt 
through the explicit binding energy of clusters. In the SMFM the clusters can be in excited 
states. In the lattice gas model also the clusters can be in excited states (excited states 
are not necessarily included in the same fashion in both models but provisions are kept for 
their inclusion). Surface tension is included in the SMFM in the binding energy relation for 
clusters (see below). This appears naturally in the lattice gas model as nucleons at the edge 
of a cluster will see an anisotropic neighbour distribution. It is because of interactions that 
the multiplicity changes in the SMFM when either the temperature or the density varies. 
The potential energy will also change with density. We therefore postulate that we have 
an interacting system where interactions lead to clusters with associated binding energy. 
Because of the short range of the interactions two different clusters do not interact except 
through the excluded volume effect. We can then apply the law of partial pressures to obtain 
an equation of state. 



IV. THE EQUATION OF STATE IN THE SMFM 

For a given freeze-out volume and temperature we can have rii nucleons, n2 clusters of 
size 2, 77-3 clusters of size 3, etc. . . An extraor dinar ly large number of choices are possible and 
the most probable distribution is obtained for that choice of ni, ?t,2, "",3 ... such that the free 
energy is minimized. We take this most probable distribution as the appropriate distribution 
for the calculation of pressure. Following the work of Ref. 0, the number of clusters of size 
k in volume Vf is given by 

nfc = exp [kPfx] ^(27rmr)^/2p/2 [pF{k)] , (4.1) 

where m=mass of a nucleon, -F(l) = and for k > 1 

F{k) = {Wk - a{T)k^'^)- <E* > +TS (4.2) 

where W = 16 MeV is the volume energy term, a{T)k'^^^ gives the surface tension correction, 
< E* >= n'^T^ / (iep) is the average excitation energy for the composite and S, the entropy, 
is S* = Ti'^T / {2ep). For the numerical work, we use (4/7r^)eir = 16 MeV. The chemical 
potential /x in Eq. (|4.1| ) is fixed from the condition ^^fc = n. The pressure is calculated 
from P{p,T) = T/Vf{J2^k)- Here p = n/{Vf + Vex)- In the above we have used the grand 
canonical ensemble. The surface tension cr(T) is taken as a function of temperature [^. For 
the discussion to follow the surface tension plays a crucial role. The SMFM-as the lattice gas 
type models-is capable of producing a U shaped distribution for Y{Z), the fragment yield 
function, against Z, the fragment charge. This would not be possible without the surface 
tension term. We also studied a SMFM model in which surface tension is not a function of 
temperature. The major features do not change in what we desribe below. 
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The P — p diagrams using this formahsm were done for A = 85 P], and also for A = 137 
(and Z = 57), a fragment size roughly corresponding to the projectile fragmentation of a gold 
nucleus. In Figs. |I] and we have shown the case of A = 137. The P — p diagram of Fig. |1| 
is strongly suggestive of a liquid- gas type phase transition. Isotherms at T = 5 and 6 MeV 
show regions where pressure is very fiat, over large variations of p. A totally constant value 
of P against p in some region implies that isothermal compressibility k = 1/ p{dp/dP)T 
goes to infinity as will happen if one moves along the Maxwell construction for a Van 
der Waals gas. This fiat region disappears above the critical point. Above the critical 
point the compressibility is a monotonically decreasing function of density. In Fig. |T], the 
pressure against p is not strictly constant in the mixed region but has a slight rise. Thus 
the compressibility is not infinite but reaches a maximum. This is sensible: infinities wil get 
replaced by maxima in a finite system. We identify the critical point by the disappearance 
of the maximum in k, (Fig. |^). For a system of mass 137 we estimate = 6.8 MeV and 
Pc = 0.4po- For a system of 85 particles the same criteria gave Tc = 6.2 MeV and pc = 0.14po- 
Values of critical temperature and density for mass number 137 are not a great deal different 
from what is obtained in the lattice gas model for an infinite system. A binding energy of 
16 MeV in nuclear matter requires a value e = —5.33 MeV. The parameters of the critical 
point in the lattice gas model will then be Tc = 6.0 MeV and pc = Po/2. 



V. THE EQUATION OF STATE IN THE LATTICE GAS MODEL 

We would like to draw P — V diagrams similar to above using the lattice gas model. 
Unfortunately, this is prohibitively difficult. While energy, fluctuations in energy (essentially 
Cv) and several other quantities can be obtained using Monte-Carlo simulations (later we 
will show results of such calculations) calculation of pressure requires evaluation of the 
partition function which can not be done using Monte-Carlo. We can draw a P — V diagram 
using mean field theory (this was already done in Refs. [1,0) but this has in some part of 
the diagram a wrong slope of pressure against p (in fact one uses that as a signature of 
phase-transition) and does not compare well with fig. 1. 

While better calculations for pressure are possible, for a limited comparison we do the 
following calculation. We use a lattice of 3'^ (this approach can not be improved by a 
straightforward increase of lattice dimensions: computation with a 4^ lattice will take 2^^ 
times as much time and is clearly beyond our scope) but use periodic boundary conditions 
to minimize edge and finite particle number effects [Q. This represents an improvement 
over the Bethe-Peierls calculation done in Ref. [^, as our lattice is in fact much larger 
than a Bethe-Peierls block. We use PV = TlnZgr, where the grand partition function is 
Zgr = e^'^Z{n). The value of A is chosen to give the required value of p/po =< n> /27 
where < n >= dlnZgr/dX. Finally V = 27a^ where = 1.0/po. In figures 3 and 4 we have 
plotted p — V and k, — p diagram. Similarities with Figs. |I] and ^ are obvious. A critical 
temperature between 5 and 6 MeV is suggested but one should keep in mind the limitations 
induced by the smallness of the lattice. The value of e used for these diagrams was —5 MeV. 
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VI. SPECIFIC HEAT 



The simple SMFM model as used here gives a maximum in as one moves across Tc 
along a density close to the critical density. This is shown in Fig. |^. This of course is in 
keeping with the findings of previous detailed investigations [|7|.|13|. Similar maxima in Cy 



are found also in the lattice gas model. We will show our results for n=137 which roughly 
corresponds to the projectile fragmentation of gold in the EOS TPC experiment. The lattice 
sizes are taken to be iV = 7^ and iV = 8^ which would correspond to two different freeze- 
out densities. For this calculation we use Monte-Carlo sampling based upon the Metropolis 
algorithm. At a given temperature we generate an event in which the potential energy is 
obtained from the number of nucleon-nucleon bonds and where the kinetic energy is simply 
the sum of kinetic energies of the nucleons. The average of many such events gives the 
average energy per nucleon, e, at a given temperature. One can then obtain per nucleon 
by numerically differentiating de/dT or else, since the events are generated in the canonical 
ensemble in the Metropolis algorithm, one can use 

C^ = ±-^[<E^>-<E>') . (6.1) 

Note that < E >= en. 

We use this opportunity to introduce an improvement to the lattice gas model. In our 
previous work we used the same value of e for attraction between like (nn and pp) and unlike 
(np) particles. Although this already gives a rather reasonable overall fit to Y{Z) against 
Z one gets a theory with di-neutrons di-protons etc. . . This can be avoided by postulating 
two kinds of bonds: that between the proton and the neutron which we denote as e„p and is 
attractive, and that between identical particles (neutron-neutron and proton-proton) which 
we denote by e^n and which is either zero or repulsive. In accordance with a molecular 
dynamics potential [|l4l which we will refer to in the next section we have used a slightly 
repulsive interaction for e„„ =1 MeV. The attractive interaction is set at e„p = — 5MeV 

In each Monte-Carlo simulation we also obtain the cluster decomposition. After a large 
number of simulation one obtains Y{Z) against Z. We expect that near the critical point 
a power law will emerge : Y{Z) oc Z~'^. The exponent r will be minimum at the critical 
point. We deduce an effective value of r from our simulation even when far from the critical 
point by using the formula ||15||: 



El'^ZYiZ) T.fzz 



— T 



— T 



(6.2) 



The variation of r against the temperature is shown in Fig. ^. The calculated against 
temperature is shown in Fig. ^ We find that the maximum in the value of the specific 
heat happens at the same temperature where the value of r minimizes. In so far as the Cy 
maximizes at the critical temperature, the prediction of the lattice gas model is similar to 
that of SMFM. We also note that using different bonds between like and unlike particles 
(which is more realistic than the simpler case) actually enhances the peak. 

In detail the calculated values of Cv's differ, specially at low temperature. Here the 
SMFM is clearly more realistic. At low temperature, the breakup event will have one large 
cluster and the specific heat of this cluster will dominate the net contribution. In SMFM the 
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low energy excitation energy characteristics is put in by hand and hence Cv has been forced 
to go hke T. In the lattice gas model the low temperature value of the specific heat per 
particle will simply be 3/2. This comes from the classical simulation of the kinetic energy 
of the constituent nucleons. Stated another way, the excited states of each cluster is put in 
differently in the two models. 

One can use the lattice gas model for predicting a caloric curve (temperature plotted 
against excitation energy), provided one assumes that calculations using the lattice gas 
model are valid beyond a certain temperature only. The uncertainty below this temperature 
can be absorbed in one constant. In Fig. [7| we show a caloric curve obtained in the lattice 
gas model. In the figure we also show a comparison with experimental data. To produce 
the caloric curve we use 

e(T) = /* C^dT' + r C^dT' = e+ r C^dT' (6.3) 
Jo Jt Jt 

We ignore the prediction of the lattice gas model below a temperature t (taken in this work to 
be 3.2 MeV) and we thus introduce a parameter e. The caloric curve beyond the temperature 
3.2 MeV is shown in Fig. ^ In the ALADIN results |T6|, T against the excitation energy 
is considerably flatter at temperature around 5 MeV, implying a large increase in the value 
of the specific heat at this temperature. We do not get such a large value of in our 
calculations. The appearance of a possible plateau in the caloric curve is also less dramatic 
in the recent EOS TPC data [p!7| . 

We end this section by summarizing the comparison between the lattice gas model and 
SMFM. As we have shown there are many common features in the two models. The lattice 
gas model has the advantage of having considerable recognition in statistical physics. It 
defines a two body interaction which is used for all its predictions. Phase transitions in the 
model are well-studied and well-defined. Another advantage is that fragment formation in 
non-spherical geometries can be studied. Such geometries do occur in heavy-ion collisions 
||T8| . The advantage of SMFM type models is that binding energy, surface tension and excited 
states can all be inputs in the model. This makes accurate comparison with experiments 
more feasible. No detailed theory of phase transitions in this type of models has yet been 
completed. It is not clear how non-spherical geometries might influence the results. 

For completeness, we now turn to a discussion of how classical molecular dynamics 
predictions differ from those of the previous approaches. 



VII. NUCLEAR PROPERTIES IN CLASSICAL MOLECULAR DYNAMICS 

Here, we choose to study equilibrium properties of nuclear systems using the methods of 
classical molecular dynamics. Our first task is to specify the nucleon-nucleon interaction. We 
use a neutron-proton interaction, fnp(r), which is attractive for large values of the separation 
distance r, but repulsive for small values. The like-particle nuclear potentials, Vnn and Vpp 
are taken to be identical and purely repulsive. Those features are chosen to satisfy the basic 
requirements of nuclear phenomenology. Our interparticle nuclear potentials are given as a 
combination of Yukawa interactions |]T1| : 

Vnn{r<rc) = Vo [exp(-/ior)/r - exp(-/iorc)/rc] , (7.1) 
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i;np(r < rj = Vr [exp(-/irr)/r - exp(-/irrc)/rc] 

- Va [exp(-/iar)/r - exp(-/iarc)/rj • (7.2) 

In the above, is fixed at 5.4 fm, a simple and adequate cutoff. It is also implied that 

Vnn{r > = Vpp{r > Tc) = . (7.3) 

The parameters of the potentials are vq = 373.118 MeV fm, Vr = 3088.118 MeV fm, Va 
= 2666.647 MeV fm, /io = 1.5 fm~^, fi^ = 1.7468 fm~^, and fia = 1-6 fm~^. The nuclear 
compressibility coefficient, K, is 250 MeV for this choice of parameters 0]. In the case of 
proton-proton interactions, the Coulomb interaction can be added separately. However, in 
order to make a consistent comparison, the Coulomb effects were left out in this work. 

In the molecular dynamics simulation, the particles are propagated in phase space by 
integrating Newton's equations of motion through a "leap-frog" algorithm [jl9|. Because we 



shall address shortly the issue of specific heat, we chose to perform simulation of canonical 



ensembles. One can show [|T9| that constant temperature dynamics are generated by the 
following equations of motion for the position and the momentum pj of each individual 
particle i of mass m subject to a force f^: 

Ti = pi/m (7.4) 
Pi = fi - (7.5) 

The quantity ^ acts as a "friction coefficient" that is recalculated at each time step and 



keeps the temperature constant to a high level of precision |20]. The system is initialized in 
momentum space from a Monte Carlo sampling of a Boltzmann distribution, and then cooled 
to the desired temperature. This method ensures a proper thermodynamic equilibrium 
distribution for the nucleons in interaction. The nucleons are confined to box whose volume 
is adjusted to yield the desired density. 

In our canonical ensemble simulations, the pressure can then be calculated for a system 



of n interacting nucleons by means of a virial expansion ||2T 



Note that we are using units where the Boltzmann constant, ks, is unity. The potential 
energy between two nucleons i and j, separated by a distance Vij, is written as Vij. It is 
straightforward to calculate the average total energy and the average total energy squared, so 
that the specific heat, C^, can easily inferred from Eq. (|6.1| ). Note that microcanical molec- 
ular dynamics has been used in the recent past to study the multifragmentation behaviour 



of heavy ion systems [22,15 



On Fig. ^ we plot pressure isotherms obtained in molecular dynamics of our nuclear 



fluid with the potential of Eq. (|7.2|) , for a system of A = 137 and Z = 57. We restrict 
our calculation to the range of temperatures and densities related to this comparative study 
with the lattice gas model and SMFM. In this temperature range, one realizes that the 
molecular dynamics predictions are devoid of the plateau regions characteristic of the above 
two approaches. It appears that the critical behaviour of the molecular dynamics ensembles 
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is in fact at a temperature just below 1 MeV and at nuclear densities around ~ po/S. Those 
parameters follow from our own simulations and are also in line with those found in a 
similar effort, using however a slightly different potential |T^. Turning to the calculations 
of Cv plotted on Fig 0, we come to conclusions similar to the ones reached in our study of 
the equation of state, namely that the specific heat does not exhibit any striking behaviour. 
There are no noticeable peaks that would follow from a flattening of the caloric curve, at 
the temperatures and densities dictated by our comparative study. At high temperatures, 
Cv tends to its value for noninteracting gases: = 3/2. We will not insist further here on a 
detailed molecular dynamics study of the thermodynamic properties of the nuclear fluid at 
finite temperatures, as this does not pertain directly to the work at hand. This is however 



done and will appear elsewhere 20 



VIII. CONCLUSIONS 

We have shown that some popular theoretical approaches to nuclear multifragmentation 
shared several common features. Those models are the percolation, lattice gas and statistical 
multifragmentation. The connection between the lattice gas and percolation models is more 
apparent and has been pointed out before, but two models that previously appeared rather 
unrelated, the lattice gas and SMFM, have some remarkable similarities to each other. We 
see this point as being quite important, as it represents a rather satisfying unity aspect from 
a theoretical perspective. This is the aspect we chose to insist on in this work. 

We also have shown some comparisons with calculations done in a classical molecular 
dynamics scenario. The results there are quite different. There are of course obvious defi- 
ciencies associated with the use of a classical approach to a system in a regime where the 
quantum aspects are undoubtedly important. We plan to investigate those deficiencies, and 
their effects, further. 
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FIGURES 



FIG. 1. A plot of pressure against p/po for the SMFM used in this work. The flat parts for the 
5 and 6 MeV isotherms signify mixed phases. 

FIG. 2. The compressibility k = 1/ p{dp/dP) at constant temperature for temperatures in the 
vicinity of the critical temperature. The maximum disappears at Tc- 

FIG. 3. The P — p diagram from the lattice gas model obtained by using a 3'^ lattice and 
imposing periodic boundary conditions. 



FIG. 4. The isothermal compressibility k at various temperatures for the 3'^ lattice. 

FIG. 5. The specific heat per particle for the SMFM used in the text for a system of 137 
nucleons. Two different freeze-out densities are used. The specific heat is dimensionless here as 
the unit for both heat and temperature is MeV. 

FIG. 6. The exponent r calculated in the lattice gas model for a system of 137 particles. We 
show results when the same value of e is used for like and unlike particles and also when they are 
different. In the latter case we show results for two densities. 

FIG. 7. The specific heat per particle in the lattice gas model for a system of 137 particles. 
The maxima here closely coincide with the positions of the minima in Fig. 6. 

FIG. 8. The caloric curve calculated according to lattice gas model for a system of 137 particles 



and compared with ALADIN data (data set II), Ref. [16| and with EOS TPC data (data set I), 



Ref. PI 



FIG. 9. Pressure isotherms calculated in classical molecular dynamics, for temperatures ranging 
from 1 to 4 MeV. We limit ourselves to densities lower than the equilibrium nuclear matter density, 
Pq. The error bars are not shown. The errors are negligible below p/po ~ 0.5, and slowly rise to 
be w 5% near p/pQ = 1. 

FIG. 10. The specific heat per nucleon, Cv, calculated in classical molecular dynamics. Three 
representative values of nuclear density are shown. The errors are statistical. 
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